Comparison between particle and fluid approximations to dust dynamics 
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Abstract: We present a new particle module of the magnetohydrodynamic (MHD) Piernik code. The original multi-fluid grid code 
based on the Relaxing Total Variation Diminishing (RTVD) scheme has been extended by addition of dust described within the particle 
approximation. The dust is now described as a system of interacting particles. The particles can interact with gas, which is described as a 
fluid. In this poster we introduce the scheme used to solve equations of motion for the particles and present the first results coming from 
the module. The results of test problems are also compared with the results coming from fluid simulations made with Piernik-MHD code. 
The comparison shows the most important differences between fluid and particle approximations used to describe dynamical evolution of 
dust under astrophysical conditions. 
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1. Piernik-MHD 

Piernik is a multi-fluid grid MHD code based on the RTVD 
conservative scheme by |Jin & Xin |1995) and |Trac & Pen 
(20031. Piernik can be used to examine dynamics of ion- 
ized or neutral gas, as well as dust treated as a pressureless 
fluid. The code computes conservative fluid variables (fluid 
density, momentum, total energy density) for each cell of 
the grid. The basic scheme has been extended by addition 
of many facilities which are useful in astrophysical fluid- 
dynamical simulations, e.g. shearing-box boundary con- 
ditions, Ohmic resistivity module and selfgravity module. 
See |Hanasz eTaLl^2008a|blc|[2009) for more details. 

2. Particle module 

Dust can be described in fluid and particle approximations. 
In the particle module of PIERNIK-MHD code the dust 
component is described as a system of independent parti- 
cles that can interact with each other. The particles can also 
interact with gas considered as a fluid. For each particle, 
equation of motion is solved using the scheme described in 
the next subsection. 

Scheme 

To solve the equation of motion for dust particles we use 
the scheme known as Verlet leap-frog method. In this algo- 
rithm, the velocities are calculated at time t+idt and used 
to calculate the positions, r, at time t + dt. In this way, the 
velocities leap over the positions, then the positions leap 
over the velocities. 

Generally, the scheme can be noted as: 



The first test problem relies on the analysis of ID sinu- 
soidal velocity perturbation. The fluid approximation re- 
sult, which can be veryfied by an analytical solution of the 
Burger's equation |Toro[ 19991, displays a conversion of the 
initial sinusoidal velocity profile into the sawtooth profile 
and then smoothing until the flat profile (figure [TJ. The dis- 
continuity in the velocity profile can be interpreted as shock 
front. 




Fig. 1 - 



Sinusoidal wave simulation result coming from the 
fluid approximation 



In the case of noninteracting particles the particle model 
leads to multiple velocity values in the velocity profile (fig- 
ure]^. To avoid the unphysical evolution of the particle sys- 
tem we have introduced interaction between particles. The 
interaction is analogous to inelastic collisions. The particles 
stick when they meet each other in the same grid cell. 




r(t + dt) = r(t) + v(t + -dt)dt, 



(1) 



Fig. 2 - Sinusoidal wave simulation result coming from the 
particle model in case of noninteracting particles 
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v(t + -dt) = v(t - -dt) + a(t)dt. (2) 

3. Results 

To compare fluid and particle approximations applied for 
the dust component we carried out several test problem 
simulations with the same initial conditions applied in both 
approaches. 
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The result coming from the particle approach with inter- 
actions taken into account (figure [3} appears different than 
the result given by the fluid approximation, because the 
particles group together into clusters. In the fluid simula- 
tions all the physical quantities are computed for every cell 
of the domain, even if density is very small. In the particle 
simulations the values of physical quantities are specified 
only in the particles locations. The fluid density profile at 
the end of the fluid simulation is represented by one peak of 
density. Respectively, at the end of the particle simulation 
all the particles are grouped together into one aggregate. 
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Fig. 3 - 



Sinusoidal wave simulation result coming from the 
particle model 



Another comparison test of the two methods is related 
to the evolution of the ID square velocity perturbation. The 
results are displayed in the figure [4] We note that the bor- 
ders of rarefying shock front area are different in the both 
metodhs. The profile given by the fluid simulation is not 
consistent with the analytical solution of the Burger's equa- 
tion. It turnes out that the difference in the results is a nu- 
merical artefact related to discontinuities in the first deriva- 



tive of the so called freezing speed function (see Trac & Pen 
(20031. Freezing speed is a quantity specific for the RTVD 
method, used to decompose vectors of conservative vari- 
ables and their fluxes into left-moving and right-moving 
waves. Freezing speed can be computed locally as veloc- 
ity for each cell or globally as the maximum velocity in 
the whole domain. When we introduced the local freezing 
speed smoothing (or just used the global freezing speed), 
the dust velocity profile became significantly more similar 
to the velocity profile obtained in the particle approxima- 
tion (the lower plots in figure [4J. 



ik-MHD with local freezing speed no smoothing 



-article module 




Piernik-MHD with global Ireezing soeed 



then in MHD with local Iteezing soeed 




Fig. 4 - Comparison of square wave simulations results 

The differences between the fluid and particle ap- 
proaches are much more significant in 2D tests. We present 
results of two identical dust fronts collision test. The left 
front moves initially two times faster than the right one. 
The fronts fragmentation (figure [5j is an effect of the parti- 
cle simulation. The fragmentation is probably caused by the 
finite grid resolution and the fluctuations of the mean mo- 
mentum in the cells containing a small number of particles. 

In the fluid simulation (figure |6| the fronts merge and 
move together with the velocity resulting from the momen- 
tum conservation law. On a longer timescale the front dif- 
fuses over the whole computational domain. 
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Fig. 5 - 



Dust fronts collision simulation result coming from 
the particle module 




Fig. 6 - Dust fronts collision simulation result coming from 
the fluid approximation 
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